Nonlinear clustering in models with primordial non-Gaussianity: 

the halo model approach 
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We develop the halo model of large-scale structure as an accurate tool for probing primordial 
non-Gaussianity (PNG). In this study we focus on understanding the matter clustering at several 
redshifts in the context of PNG that is a quadratic correction to the local Gaussian potential, 
characterized by the parameter /nl- In our formulation of the halo model we pay special attention 
to the effect of halo exclusion, and show that this can potentially solve the long standing problem of 
excess power on large scales in this model. The halo model depends on the mass function, clustering 
of halo centers and the density profiles. We test these ingredients using a large ensemble of high- 
resolution Gaussian and non-Gaussian numerical simulations, covering /nl = {0, +100, —100}. In 
particular, we provide a first exploration of how halo density profiles change in the presence of PNG. 
We find that for /nl positive/negative high mass haloes have an increased/decreased core density, 
so being more/less concentrated than in the Gaussian case. We also examine the halo bias and show 
that, if the halo model is correct, then there is a small asymmetry in the scale-dependence of the bias 
on very large scales, which arises because the Gaussian bias must be renormalized. We show that 
the matter power spectrum is modified by ~ 2.5% and ~ 3.5% on scales k ~ 1.0 ftMpc" 1 at z = 0.0 
and 2 = 1.0, respectively. Our halo model calculation reproduces the absolute amplitude to within 
< 10% and the ratio of non-Gaussian to Gaussian spectra to within < 1%. We also measure the 
matter correlation function and find similarly good levels of agreement between the halo model and 
the data. We anticipate that this modeling will be useful for constraining /nl from measurements 
of the shear correlation function in future weak lensing surveys such as Euclid. 
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I. INTRODUCTION 



Over the last few decades, through experiments such 
as 2dFGRS, SDSS and WMAP [l-i], great strides 
have been made in quantifying the parameters of 
the perturbed and unperturbed Friedmann-Lemaitrc- 
Robertson- Walker (FLRW) model. Besides the detailed 
information about the present Universe, these experi- 
ments have also opened up important new windows into 
the physics of the early Universe. 

The inflationary paradigm is so far the leading physi- 
cal explanation for the origins of structure. The single- 
field slow-roll theory makes four fundamental predictions: 
a flat universe (quantified by the curvature density pa- 
rameter fix); a primordial density power spectrum with 
power-law index (n s ) just less than unity; a nearly Gaus- 
sian distribution of primordial density fluctuations (de- 
viation from Gaussianity being quantified by the /nl pa- 
rameter); a spectrum of gravitational waves (character- 
ized by the amplitude of the quadrupole tensor to scalar 
ratio r). The measurement of the CMB temperature 
anisotropics combined with a number of cosmological 
probes, such as the Baryonic Acoustic Oscillation (BAO) 
scale in the SDSS LRGs, and the present day value of 
the Hubble constant, provide strong supporting evidence 
in favor of the first three of these predictions. The cur- 
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rent WMAP7+BAO+i?o combined constraints on these 
parameters are 0: -0.0133 < n K < 0.0084 (95% CL); 
n s = 0.963 ± 0.012 (68% CL); / NL = 32 ± 21 (68% CL). 
For the spectrum of gravitational waves current CMB 
experiments (WMAP7+ACBAR+QUaD) place an up- 
per bound of r < 0.33 (95% CL) (j-i]. The latter can 
only be further falsified with dedicated CMB polarization 
experiments such as PLANCK 0. 

Intriguingly, the constraints on the amount of primor- 
dial non-Gaussianity (hereafter, PNG) from the CMB arc 
tightening around a non-zero value, there currently being 
~ 1.5<7 evidence against pure Gaussian scalar perturba- 
tions [1] (see also @ who claimed a 99.5% CL detection of 
non-zero /nl)- If /nl is found to be substantially greater 
than zero at high significance, then this would rule out all 
inflation models based on a single scalar field [9( . On the 
other hand, multi-field models may produce large levels 
of PNG and also scale-dependent /nl- Thus testing the 
Gaussianity of the initial fluctuations is of prime concern. 

Most tests for PNG have primarily been performed 
on the temperature anisotropics in the CMB, however 
for several decades it has been theoretically understood 
that PNG also affects a number of large-scale structure 
observables However, very few tests with real 

galaxy survey data have been performed. This was pri- 
marily due to the fact that in the present day Universe 
the density fluctuations do not remain pristine, but have 
been driven to a non-Gaussian state by gravity. Gravi- 
tational evolution of the density perturbations correlates 
the amplitudes and phases of different Fourier modes, 
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thus one is faced with the problem of decoupling primor- 
dial from gravitational non-Gaussianity. Furthermore, in 
observing large-scale structures one does not in general 
get information about all points in space, but instead one 
is restricted to learning only about the galaxy distribu- 
tion. In general this is related to the underlying density 
statistics through a bias function (or functional), which 
may be complicated and stochastic. In the simple case 
of deterministic local biasing, one may attempt to use 
higher order statistics such as the galaxy bispectrum to 
disentangle the effects of gravity and galaxy bias [bil - [l7j . 
However, no constraints on /nl have yet emerged from 
such schemes. This partly owes to the fact that past 
survey volumes have been too small for such tests to be 
performed with any confidence. The survey volumes of 
ongoing and future planned missions will surely change 
this. 

Recently, in a ground breaking paper, it was shown 
both theoretically and in numerical simulations by [l8j . 
that there is a strong signature of /nl in the power spec- 
trum of dark matter haloes. The effect is to induce a 
scale-dependent bias correction Ab cx k~ 2 . The excit- 
ing prospect of this is that, since the signature affects 
primarily only the largest scales, k < 0.02 /iMpc -1 , one 
can in principle decouple the effects induced by nonlinear 
bias and gravity from those of PNG, and so constrain 
the latter. There has been much activity in quantify- 
ing the effects of this scale-dependent bias on the power 
spectrum in simulations |19l - l2ll |: and also there has been 
much theoretical activity to truly understand how the 
scale-dependent bias arises |22H2a |. and for a current re- 
view see [27j]. This has culminated in several recent at- 
tempts to constrain /nl from large-scale structure data 

Another recent result has shown that the nonlinear 
dark matter power spectrum alone is sensitive to the pres- 
ence of / NL [SElIlillllll- This prompted [31] to pro- 
pose that the statistics of the initial conditions could be 
tested through weak gravitational lensing by large-scale 
structure, in particular through measurement of the two- 
point shear correlation functions. The advantage of such 
a probe is that it is only sensitive to the total mass dis- 
tribution projected along the line of sight. On the down 
side, the signal is weak, with deviations being of the or- 
der several percent. However, future all sky weak lensing 
surveys such as EUCLID and LSST || will be able 
to probe changes in the convergence power spectrum to 
percent level accuracy. It is therefore of great interest 
to quantify in detail how sensitive such a mission would 
be to constraining PNG of the local type. The starting 
point for such a study is an accurate model of the matter 
power spectrum as a function of redshift. 

In this study we focus on understanding how PNG 
shapes the two-point matter clustering statistics in 
Fourier and configuration space. Currently there is no 
accurate analytic model for describing the effects of PNG 
on the matter power spectrum over the range of scales 
that will be required for future weak lensing missions 



(k £ [0.01, 100.0] /iMpc" 1 [MbML Tm s we shall at- 
tempt to build. Following [3a], [31j proposed that the 
halo model might be able to predict the power spectrum. 
They suggested that PNG would modify predictions in 
two ways: through the abundance of massive clusters; 
and through the scale-dependent bias. No attempt to 
compare their model with numerical simulations was pre- 
sented and so it remains an open question as to how re- 
liable this proposition is. Moreover, they assumed that 
halo profiles are not affected by PNG. In this work we 
show, through detailed numerical simulations, that this 
assumption is wrong, and that there are important effects 
which need to be taken into account if the halo model is 
to be used to make accurate predictions. 

Furthermore, we make important new developments in 
the importance of halo exclusion in order to make robust 
predictions, and we show for the first time how this may 
also resolve a long-standing problem with the halo model, 
the excess power on large scales. 

Lastly, when exploring how the scale-dependent bias 
from PNG enters the halo model framework, we have 
shown that if the model is to be self-consistent, then 
there must be a modification to the original formula of 
[l8| . This modification creates an asymmetry in the bias 
for models with the same value of /nl but with an op- 
posite sign. This arises due to the fact that it is not the 
Gaussian bias which enters the formula of Dala, but in 
fact the total non-Gaussian scale-independet bias. We 
believe that inclusion of this effect will slightly modify 
current constraints on /nl- 

The paper is broken down as follows: In fJTT] we 
overview the local model for PNG and its impact on den- 
sity statistics. In £11111 we describe the suite of A-body 
simulations used in this study. In £HVI we explore model- 
ing of the matter power spectrum with PNG, in particu- 
lar first explore the perturbation theory approach, before 
moving on to describing the halo model. In SjV] we per- 
form a phenomenological study of the necessary ingredi- 
ents of the halo model and their dependence on PNG. In 
Q VII we present our predictions for the halo model power 
spectrum as a function of redshift and also the two point 
correlation function. Finally, in Will we summarize our 
findings and conclude. 



II. PRIMORDIAL NON-GAUSSIANITY 

A. Potential statistics 

We shall be working with the local model of PNG, 
that is we consider only local quadratic corrections to the 
Gauge invariant Bardeen's potential perturbation, which 
on scales smaller than the horizon size reduces to minus 
the Newtonian potential [la ]: 



$(x) = 0(x) + / NL [<Xx) 2 - (0 2 (x))] 



(1) 



where 4>{x) is the Gaussian potential perturbation after 
matter radiation equality, scaled in terms of units of c 2 
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to yield a dimensionless quantity. Following standard 
convention, it is defined to be related to the Newtonian 
potential as $(x) = _$Newton( x )_ The term is 

subtracted to ensure that $ is a mean zero field. In linear 
theory the typical fluctuations are of the order $ ~ 10~ 5 , 
and so the non-Gaussian corrections are of the order ~ 
O.1%(/nl/1OO)(0/1O- 5 ). 

This transformation of the Gaussian potential leads to 
a small correction to the power spectrum, but its main 
effect is to generate a primordial potential bispectrum. 
To see this consider the Fourier transform of Eq. (JT|) : 



$(k) = ^(k) + /nl 



d 3 q 

(2tt) 3 



<^(q)^(k - q) 



(2) 



The power spectrum, defined as ($(ki)$(k 2 )) 
P, I >(|ki|)(27r) 3 (5 £ '(ki + k 2 ) for this field is given by: 



P*(ki) = P^(k 1 )+2/ ] 



XL 



d 3 qr 

(2tt) 3 



P^(qi)P^(|ki-qi| 



The three point function in Fourier space is given by 

($(ki)3(ka)$(ks)) 
d 3 q 3 



(3) 



/nl 
+/nl 



(^(ki)^(k 2 )0(q 3 )^(k3 - q 3 )) + 2cyc 
(qi)^(ki - qi) 



(2tt) 3 

:i / d 3 qi d 3 q 2 d 3 q 3 



(2tt) 3 (2tt) 3 (2tt) 3 
x 0(q2)0(k 2 - q 2 )0(q 3 )0(k 3 



03, 



(4) 



Recalling that the expectation value of odd powers of 
the Gaussian variables vanish, and from Wick's the- 
orem we have that the even powers can be writ- 
ten as (0(ki) . . . 0(k„)) = £ allpairs Ut^ MkiM-kO). 
Also, defining the bispectrum as ($(ki)$(k 2 )$(k 3 )) = 
B$(ki,k 2 ,k 3 )(27r) 3 <S D (ki + k 2 + k 3 ), then we find the 
primordial potential bispectrum to be: 



B$(ki,k 2 ,k 3 



+ 4/j 



2/nl [P (k!)P (k 2 
d 3 qr 



NL 



(27T) 

^(|qi+k 2 |) 



2 eye] 
3 P*(qi)P^(|q 1 



kil 



P (|qi+k 3 |) (5) 



Restricting our attention to the case where /nl < 100 
then we may safely neglect the second terms on the right- 
hand-side of Eqs © and © [zi| 



B. Density statistics 

The primordial matter potential and density fluctu- 
ations, extrapolated to the present day, can be related 
through Poisson's equation: 

V 2 $(x,a) = -^[p(x,a)-p(a)]a 2 

3 (H \* D(a) 
2 \ c J a 



where a is the expansion factor, p(a) = Sl(a)p CI it(a) oc 
a -3 is the mean density of the Universe, fi(a) is den- 
sity parameter, p C rit(o) the critical density, H(a) is the 
Hubble parameter, D(a) is the linear growth factor nor- 
malized to be unity at the present time, and quantities 
labeled with a subscript are to be evaluated at the 
present time a,Q. This equation may be solved in Fourier 
space to give an explicit relation for the potential: 

Evolving the potential back to the initial epoch a,, and 
dividing the transfer function, then we have the following 
relation between the present day density and primordial 
potential perturbations, 



<5o(k, a ) = a(k, a % , a )$(k, a*) 



where we have defined 



a(fc,Oj,a ) 



2c 2 g(a l ,a )k 2 T(k 1 a ) 



(8) 



(9) 



In the above equation g(ai,a 2 ) = [D(a 2 )/ D(ai)] [ai/a 2 ] 
is the growth suppression factor (ai < a 2 ), and for 
LCDM, g{ai,a ) « 0.75. 

In possession of the mapping from the present day 
density to primordial potential perturbations through 
Eq. (O , we may now examine statistics of the density 
field. The most important statistic that we will need to 
know is the present day skewness filtered on the mass 
scale M. Following Appendix [Ai this can be written: 



{S 3 M (x : a Q )) = 6/nl 



2tt 



:AtfW(fci,M)a(A;i)P0(fci) 



/dk 
^k 2 2 W(k 2 ,M)a(k 2 )P^(k 2 ) 



d/j,W{k 3 ,M)a(k 3 ) 



(10) 



where k 3 



k'i + 2k 1 k 2 /j,, W(k, M) is a filter func- 



k 2 + k 2 

tion that selects the mass scale M, and where for 
brevity we shall make the definitions P^,(k 2 ) = P^(fc 2 , ai) 
and a(k) = a(k,ai,ao). In what follows we shall 
also make use of the reduced skewness, defined to be 
53 (M) = (Sl { (x)) I (S 2 j(x)) . Finally, for completeness, 
we may write the density power spectrum at some arbi- 
trary epoch a in terms of the primordial potential power 
spectrum simply as: 



Ps(k, a) = D 2 (a,a Q )a 2 (k,a il ao)P^(k 1 a i ) 
= D 2 {a)a 2 (k)P 4> {k) . 



ill) 



In Figure [T] we show the present day skewness of the 
density field S 3 (M)a 4 (M), as given by Eq. (10) and with 
/nl = 1- For the calculation we employ a real space 
spherical-top-hat filter with radius given by the mass 
scale R — (SM/Anp) 1 ^ 3 , which in Fourier space has the 
form: W(y) = (3/y) [siny — ycosy], where y = kR. The 
figure shows that the skewness is very small, < 1CP 3 for 
R ~ 5/i- x Mpc; and if / NL = 100, it is < 0.1. 
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10 13 10 14 10 15 

Mass [h _1 M ] 

FIG. 1: Statistical description of the present day non- 
Gaussian density field, with /nl = 1, as a function of La- 
grangian mass scale M. The solid line shows the standard 
deviation of the density field; the dot-dash line shows the 
skewness of the density field; the dash and dotted lines show 
the logarithmic derivatives of these quantities. 



of comoving radius R = 8 h 1 Mpc. 

The matter transfer function was generated using CAMB 
[3^ . All of the simulations were run using the A-body 
code Gadget-2 38]. The same Gaussian random seed 
field 4> is employed for each /nl = {0, +100,-100} sim- 
ulation set, and varied between sets. This allows the 
sampling variance between different models of /nl to be 
minimized when we construct statistics from the ratios 
of observables. The initial particle distribution was gen- 
erated at redshift Zi = 99 using the Zel'dovich approxi- 
mation [391 ] . 

Regarding the generation of the initial conditions for 
the non-Gaussian simulations, we adopt the standard 
(CMB) convention in which $(x) is primordial, and not 
extrapolated to present epoch. Furthermore, we point 
out that the local transformation to the potential given 
by Eq. ([1]) is performed before multiplication by the mat- 
ter transfer function T(k). 

Dark matter halo catalogs were generated for all snap- 
shots of each simulation using the Friends-of-Fricnds 
(FoF) algorithm We set the linking- length param- 
eter to the standard b — 0.2, where b is the fraction of 
the inter-particle spacing. For this we used the fast par- 
allel B-FoF code, kindly provided by V. Springel. The 
minimum number of particles for which an object was 
considered to be a bound halo, was set to 20 particles. 
With particle mass m p s» 3.0 x 10 11 h^ 1 M Q , this gave us 
a minimum host halo mass of 6 ~ 10 12 M Q //i [79] ]. 



III. THE TV-BODY SIMULATIONS 

In order to explore the impact of PNG on the clus- 
tering statistics of the density field, we have generated 
a large ensemble of high-resolution A-body simulations 
of the ACDM cosmology seeded with Gaussian and non- 
Gaussian initial conditions. This set is an augmented 
version of the ensemble used by Desjacques et al. [H} to 
study the mass function of CDM haloes and their bias, 
and in Sefusatti et al. [171 ] to explore the matter bispec- 
trum. 

The non-Gaussianity in the simulations is of the local 
form, c.f. Eq. (pj, and we use twelve sets of three simula- 
tions, each of which has /nl = 0, ±100. Each simulation 
was run with A — 1024 3 particles in boxes of side L — 
1600/i -1 Mpc and V Sim ~ 4.096 h~ 3 Gpc 3 , and this gives 
us a total simulated volume of Vrot ~ 49.152 /i^ 3 Gpc 3 . 
The interparticle forces were softened on scales of 0.04 
times the mean inter-particle distance, which corresponds 
to / so ft w 40/i _1 kpc. We used the WMAP5 cosmological 
parameters: {h = 0.7, fi m = 0.279, fl h = 0.0462, n s = 
0.96}, and a normalization of the curvature perturbations 
A^(fc) = 2.21 x 10- 9 (fc/fco)™ s "\ with k = 0.02Mpc"\ 
where the curvature perturbation is related to the scalar 
potential 1Z = 5$/3, and $ was defined in Eq. ((T|). In 
terms of the variance of matter fluctuations linearly ex- 
trapolated to the present day, this gives as « 0.81, where 
the variance is computed with a spherical-top-hat filter 



IV. MODELING THE NONLINEAR POWER 
SPECTRUM 

A. Gravitational perturbation theory 

The nonlinear evolution of the density field can be fol- 
lowed using gravitational perturbation theory. In this 
approach one writes down the equations of motion for 
the CDM fluid in an expanding spacetime and looks for 
a series expansion in the density and velocity divergence 
of the CDM fluid [for a review of the subject see|4l[. The 
important point to note is that the perturbative solutions 
do not change when we consider the case of structure for- 
mation in models with PNG. However, what does change 
is the way in which the statistics of the density field be- 
have. This owes to the fact that there is now a complete 
hierarchy of connected correlators of the field. Following 
the work of (24[, and keeping only terms that are lin- 
ear in /nl> one finds that up to the next-to- leading-order 
corrections to the matter power spectrum can be written: 

P(k) = P n (k) + P 22 (k) + P 13 (k) + P 12 (k) , (12) 

where Pn, is the linear power spectrum, and Pi 3 and 
P22 are the so-called one-loop corrections, which appear 
in the standard Gaussian theory [4l| . The new term for 
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FIG. 2: Comparison of the nonlinear matter power spectra measured from the suite of ./V-body simulations with predicitons 
from nonlinear gravitational perturbation theory. The left and right panels show results for z — 1.0 and z = 0.0, respectively. 
The top section of each panel presents the absolute power, the bottom section shows the ratio of the nonlinear power in the 
non-Gaussianity to the Gaussian models. Points with errors denote estimates from the simulations and lines denote theoretical 
predictions. The colors green, red and blue denote the models /nl = {0, 100, —100}, respectively. 



PNG is P12, which has the form [Hj]: 



non-Gaussian and Gaussian models as: 



2£ 3 (a)/ N La(fc)P (fc)fc 3 



7(2tt) 

d y .3, 



-^-y i a{yk)P 4 ,{yk) 
>y mi n V 

x / d^a{kip) 

J max[— l,ei] 



1 + V 2 - 2y/i 



1 



P^{k) P^yk) 



where tfj 2 = k 2 (l + y 2 — 2y/z), y n 



/ k and y„ 



(13) 



k m in/k. We also defined t\ = [k 2 + q 2 — k 2 lax \ /2kq and 
e 2 = [k 2 + q 2 — fcmin] l^-kq. The cut-off scales are set to 
be /cmax = 10/iMpc" 1 and fc m ; n = 2tt/L, with L being 
the simulation box length. Note that the spectrum P\i 
arises due to the existence of the non-zero primordial 
potential bispectrum (c.f. Eqn. [5j. 

Let us now define the ratio of the power spectra in the 



/3pT(fc,/ NL ) = 1 + 



Pn(k,a) + P 13 (k,a) + P 22 (k,a) 



(14) 

In [Tg it was shown that, for k < 0.2/iMpc -1 , the PT 
description was able to capture to high precision the same 
ratio measured from iV-body simulations. In Figure[2]we 
extend the analysis to much higher wavenumbers, rele- 
vant for cosmological weak lensing studies. On the left- 
hand-side of the figure, we show the results for z = 1.0 
and on the right those for z = 0.0 (we provide further 
details of these measurements in WI|) . The top sections 
of the panels show the absolute power and on a log-log 
plot the points for the three /nl models cannot be dis- 
tinguished. The bottom sections show /3pt(A;, /nl)- We 
see that the effect of /nl = ±100 is to induce ± ~ 3.5% 
modulations in the nonlinear matter power spectrum at 
z = 1.0 and on scales of order k = 1.0/iMpc , and that 
this reduces to ~ 2.5% by z = 0.0. These deviations, 
although small would be larger than the measurement 
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errors in future weak lensing missions and so need to be 
accurately characterized. 

The figure also shows that, whilst PT captures well 
the behavior of /3px measured in the simulations on very 
large scales k < 0.2/iMpc -1 , it fails to model the re- 
sults on smaller scales. It is also worth pointing out that 
whilst /3px is well characterized for k < 0.2 h Mpc -1 , 
the absolute power does not match well the measured 
absolute power at that scale, it being a factor of ~ 2 
higher than the measurements at z = 0.0. This leads us 
to explore alternative approaches to modeling the non- 
linear power spectrum on smaller scales. Before mov- 
ing on, we note that there are interesting new analytic 
approaches to solving the pertubation theory at higher 
orders being developed, i.e. renormalized perturbation 
theory and renormalization group theory |42h44| . These 
will push our analytic understanding deeper into the non- 
linear regime k > 0.2 h Mpc" 1 . However, we point out 
that the most recent calculations in the literature, show 
only a very modest improvement over the standard PT 
that we considered in this paper 143;]. Moreover, whilst 
certainly interesting, these techniques are unlikely to be 
able to probe the range of wavenumbers required for fu- 
ture weak lensing missions (k € [0.01, 100.0] /iMpc -1 ). 



B. The halo model approach 

The halo model was developed by a number of authors 
[4I5I - I47I ] . and for a review see [HJ]. In this model all of 
the mass in the Universe is distributed into dark matter 
haloes, each labeled with some physical properties. Typ- 
ically one simply labels each halo by its mass, however 
more complicated approaches can take into account, for 
instance, halo shape [491 ] . 

In the halo model the density field of dark matter may 
be written as a sum over all haloes, 



N 



/9 (x)=^M i C/ i (x- 



\M 3 ) 



(15) 



where N is the total number of haloes, Mj and x^- are 
the mass and center of mass of the jth halo and Uj = 
Pj(xj)/Mj is the mass normalized density profile. The 
statistics of the density field may be computed directly. 
In particular for the density power spectrum it can be 
shown that it can be written as the sum of two terms [for 
more details see W3, KM : 



particles in distinct haloes. They have the explicit forms: 
1 f 00 

Pm(k) = / dMn(M)M 2 \U(k\M)\ 2 ; (17) 
P Jo 

-1 />oo 2 

P2H(k) = — / T\{dMin(Mi)MiUi(k\Mi)} 
P Jo ;=1 

x P^ nt (k\M u M 2 ) , (18) 

where the essential new ingredient is the power spec- 
trum of halo centers with masses Mi and M 2 , denoted 
P^ t (k\M u M 2 ). 

C. The halo-center power spectrum 

The power spectrum of halo centers contains all of 
the information for the inter-clustering of haloes; pre- 
cise knowledge of this term is required to make accurate 
predictions on both large and small scales. In principle, 
P^ t (k\Mi, M 2 ) is a complicated scale-dependent func- 
tion of Mi, M2 and k [5lJ. The usual starting point for 
modeling this is to assume the local deterministic bias- 
ing approach (52l.[53j. This can be summarized as follows. 
Consider the density field of all haloes with masses in the 
range M to M + dM, smoothed with some filter of scale 
R. We now assume that this field can be related to the 
underlying dark matter field, smoothed with the same 
filter, through some deterministic mapping and that this 
mapping should apply independently of the precise posi- 
tion, x, in the field: 



Q (x|ii,A/)= F {M . R} ]<5(x|i?)] 



(19) 



where the sub-scripts on the function T indicate that it 
depends on the mass of the haloes considered and the 
chosen filter scale. The filtered density field is 



6(x\R) = ^J d 3 y 5(y)W(\x-y\,R) 



(20) 



W(|a;|,i?) being some normalized filter. Taylor expand- 
ing J-{m,r} about the point 5 — yields 



T {M , R} [6(x\R)] = 



h(M\R) 



[S( X \R)Y 



(21) 



P(k) = P m (k) + P 2H (fc) , 



(16) 



We now assume that there is a certain filter scale above 
which J-{m,R) is independent of both the scale considered 
and also the exact shape of the filter function. Hence, 

8\ x \R,M) = ^ h -^[5{K\R)Y , (22) 

i=0 



where the first term is referred to as the '1-Halo' term, 
which describes the intra-clustering of dark matter par- 
ticles within single haloes; the second term is referred 
to as the '2-Halo' term, and describes the clustering of 



where the bias coefficients are 



dx l 



(23) 
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The bias coefficients from the Taylor series are not 
independent, but obey the constraint, (£ h (x|i?)) = 0, 
which leads to 



bn {§n) 



(24) 



Thus, in general bo is non-vanishing and depends on the 
hierarchy of moments. This allows us to re- write equation 
as 

5\ X \R,M) = £ b -^fl {[6(x\R)f - <*(*!*)>} • 

l — l 

(25) 

Nevertheless, we may remove &o from further consider- 
ation by transforming to the Fourier domain, where it 
only contributes to 5(k = 0). 

Thus the halo centre correlation function has the form: 

& r ) = (6 h (x\R,M)6 h ( X + r\R,M)) 

= b 1 (M 1 )b 1 (M 2 )^ R (r) + I [b 1 (M 1 )b 3 (M 2 ) 

6 

+ 63(Mi)6i(M 2 )] (J«(x)4(x + r)) 

b 2 {M 1 )b 2 (M 2 ) , a 2 
+ 1 V 5 R( x ) 5 fl( x + *)) + ••• (26) 

where £r(t") is the nonlinear matter correlation function 
smoothed on scale R. Fourier transforming the above 
exp ression we obtain the halo center power spectrum [5ll 
MM: 



P^ nt {k\M u M 2 ) = b x (M l )b 1 {M 2 )P Nlj {k\R) + 0{b 2 , 



••), 

where the parameters bi are the nonlinear bias coefficients 
and Pysh{k\R) is the nonlinear matter power spectrum 
smoothed on scale R. It has recently been proposed that, 
for PNG, the halo bias is also a function of the local 
gravitational potential [2(1 [55[ , we shall not explore this 
possibility here, but simply note that it should give rise 
to the same scale dependence of the linear bias. Further, 
since this is a first order attempt to calculate the effects 
of PNG on the matter clustering in the halo model, we 
shall restrict our attention to the case of linear bias and so 
neglect terms 6, with i > 1. Whereupon, P^ t becomes a 
separable function of mass and scale. We present details 
of the &i(M,/ NL ) model in JV] For P NL (fe, R, / NL ), we 
make the simple approximation: 

PNL(k\R,hL) = W 2 {kR)P h8Llollt {k)p PT (k,f NL ) (28) 

where W(kR) is a smoothing function, Phaiofit(fc) is the 
nonlinear matter power spectrum model of [561 ] , valid for 
Gaussian initial conditions, and where /3px was defined 
earlier in Eq. (Tl4"|) . 

As was argued in [5l|, H3, [H[ another essential com- 
ponent of the inter-clustering of haloes is halo exclusion. 
That is, one must remove the correlations which arise on 
scales inside the sum of the virial radii of the two haloes 



M i and M 2 . As was shown in [5l| , this effect can formally 
be written: 



^ nt (r\M u M 2 ) = -1 ; (r < r vhA + r vir , 2 ) 



(29) 



where r v i r is the virial radius of a halo and where £ hh is 
the correlation function of dark matter halo centres, de- 
fined: ^ nt (r\M 1} M 2 ) = (S h (x\M 1 )S h {x + r\M 2 )). The 
— 1 in the above is simply the value that £ must obtain 
in order for the joint probability of finding halo center 
separations r < r v - lr ^ + r v - lr ^ 2 to be zero. In the litera- 
ture, various approximate schemes have been proposed to 
model the exclusion effect [13, , these involve placing 
a cut-off in the upper limit of the mass integrals in the 2- 
Halo term. We shall not follow such schemes, since these 
approaches do not reproduce the correct power spectrum 
asymptotics for the exact calculation, which we show be- 
low. Instead we follow [5l|, and evaluate the above ex- 
pression exactly. In this case the halo center power spec- 
trum can be written in terms of the correlation function 
of halo centres as: 



d a T&*(k\M u M 2 ,R)j (kr) 



(30) 



On inserting Eq. (]29[) for scales inside r v j rj i + r v i rj 2 and 
the relation ^(ftMuM^R) = 6(Mi)6(M 2 )C(t-|!r) on 
larger scales, where £,(r\R) is the dark matter correlation 
function smoothed on the scale R, then we find 



P c h c h nt (A:|i?) = / dhb(M 1 )b(M 2 )ar\R)jo(kr) 

f" vir .i+I'vii.2 

^3. 







dM-l)jo(fcr) 
d 3 r6(M 1 )6(M 2 )^(r|i?)j (A:r) 



'Vip.l+I'Vir ,2 



d 3 r [1 + b{M 1 )b(M 2 )^r\R)] j (kr) 



= P, 



NoExc.hh 



(k)-p^ : h \k) 



cent 



(31) 



The first term in the last line of the above equation rep- 
resents the usual expression for the clustering of halo 
centers, and the second term represents the correction 
due to halo exclusion: 



pNoExchh _ b ( Ml )b(M 2 )P U n(k) ; 



cent 

-,Exc.hh 



(32) 



P 



cent 



V,r,1 d 3 r [1 + 6(M!)6(M 2 )f(r)] jo(fcr)(33) 



Taking s = s{Mi,M 2 ) = r virj i + r vir;2 and y = ks, 
we may deduce the following asymptotic properties for 

pExc,hh, 
cent 

• Large-scale limit: for the case k — ¥ 0, we have that 
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lim P.. 

y^O 



jo(kr) — > 1 and so 

c J lh (k\M u M 2 ) 

= f d 3 r [1 + b(M 1 )b(M 2 )£(r)] 
Jo 

= V(s) [1 + b{M 1 )b(M 2 )l{ s )] , (34) 

where V(s) — 4irs 3 /3 is the volume of the exclusion 
sphere for the haloes Mi and M 2 and where £ is the 
volume averaged correlation function. This appears as 
white noise power contribution, and so it acts to reduce 
any large-scale shot noise component. 



Small-scale limit: for the case fc 



we have that 



lim P c e S hh (fc|Mx,M 2 ) 



lim / d 3 r[l + ^Mi)6(M 2 )f(r)]io(*r) 



fc— >oo 



= lim 



fc 3 



d 3 yjo(y) 



l 



+ Um |p ^ A z y'b{M 1 )b{M 2 )i{y l /k)j Q {y')^ 

= (2nf6 D (k)+P^- hh (k). (35) 

Thus the effect of halo exclusion on small scales, is to 
exactly null the 2-Halo term without exclusion. 

This leads us to write the full 2-Halo term as 
P 2H (k)-P 2 l c (k), where 



Pm(k) = t~ / TT {dMin(Mi)MiUi(k\Mi)} 
P 2 Jo 11 



1=1 



d A v [1 + 6(Mi)6(M 2 )€(r)] jo(Ar)(36) 



As a short aside, in Appendix [B] we forward the idea 
that halo exclusion may resolve the well known problem 
of excess large-scale power in the standard formulation 
of the halo model and that after taking this into account 
the theory is consistent with perturbation theory results, 
like RPT (Hi]. 



HALO MODEL INGREDIENTS IN 
NON-GAUSSIAN MODELS 



The model that we described in gTVBl specified noth- 
ing about the cosmological model, other than that the 
end state of gravitational clustering leads to the forma- 
tion of a distribution of haloes with some characteristic 
spectrum of masses, density profiles, and that halo cen- 
ters are clustered. Thus no extension of the formalism 



is necessary in order to use the halo model to describe 
clustering in more exotic models, such as PNG. However, 
what must necessarily change are the ingredients of the 
model: the mass function, the halo bias and the density 
profiles. We now study these in detail in the context of 
PNG. 



A. The halo mass function 

The mass function of dark matter haloes in models 
of structure formation from Gaussian initial conditions 
has been widely studied over the past decades fofT[66| . 
Conventionally, the mass function is represented: 



dn 



d\ogM 



M v 'dlogM 



Sc(z) 
a{M) 



(37) 



where for the Press-Schechter (PS) and Sheth & Tormen 
(ST) mass function we have: 



/ps(i>) 



— z^exp 

7T 



T 



A\j-y/qy[\ + (y/qv)- 2 »] exp 



2 



(38) 
(39) 



For the ST mass function the amplitude parameter A is 
determined from the constraint § d\ogvf{v) — 1: which 
leads to A- 1 = {1 + 2-*T [0.5 -p] /V [0.5]}. ST's orig- 
inal parameters are: {A — 0.3222, p = 0.3, q — 0.707}. 
The derivative term on the right-hand-side can be cal- 
culated from Eq. (|A8I) in Appendix [AJ Fig.[T] shows the 
derivative of Eq. (|A8|) computed for the WMAP5 param- 
eters compared to the skewness due to PNG. 

A number of authors have studied the effects of PNG 
of the local type on the mass function of dark matter 
haloes 0, M, M, HI HJ WH The most impor- 

tant task in extending the Press-Schechter framework 
is to find an analytic expression for the f-point prob- 
ability density function (PDF) of the smoothed matter 
fluctuations. ' 13. hereafter MVJl gave the first formal 
derivation using a path-integral approach. [60l hereafter 
LV] used the Edgeworth expansion (or more simply the 
Gram-Charlier Type la series) to recover the PDF. The 
key idea of these expansions is to write the characteristic 
function of the non-Gaussian PDF to be approximated in 
terms of the characteristic function of the Gaussian PDF, 
and to then recover the non-Gaussian PDF through the 
inverse Fourier transform method [see [69|. 

For small amounts of PNG the mass function can be 
written: 



dn{Mf Nh ) dn ST (M) 



d log M 



dlogM 



R{v, /nl) 



(40) 



where the ratio of the non-Gaussian to the Gaussian mass 
functions of MVJ and LV can be written: 
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FIG. 3: Mass functions of haloes in models with local primordial non-Gaussianity as a function of FoF (6 = 0.2) halo mass and 
at several redshifts. Left and right panels compare estimates from the ensemble of iV-body simulations with the theoretical 
predictions of Lo Verde et al. [6y| and Matarrese et al. respectively. Gaussian mass function predictions are given by 

Sheth and Tormen Top panels show the absolute mass function. The green, red and blue points with errors denote 

estimates for /nl = {0, +100, —100}. The solid, dash and dot dash lines represent the predictions as given by Eg. (|40[1 for 
/nl = {0, 100, —100}. From bottom to top the points show results for redshifts z € {1.0, 0.5, 0.28, 0.0}. Bottom panels: ratio 
of the estimated mass functions with their respective theoretical predictions. For clarity, we have offset the results for redshifts 
z e {1.0, 0.5, 0.28, 0.0} by {0.0, 0.2, 0.4, 0.6} in the positive y-direction. 




FoF Mass [h *M „] 



FoF Mass [h ^ 



FIG. 4: Ratio of the non-Gaussian to Gaussian mass functions as a function of FoF (b = 0.2) halo mass. Left and right panels 
show compare the simulations with the theoretical models of [fj^] and [l3T ]. respectively. The points with errors denote the 
ensemble averages of the mass function ratios measured at expansion factors z £ [1.0, 0.5, 0.28, 0.0], where larger point symbols 
denote later times. Solid lines denote theoretical models. The theory predictions were generated using: 5 C — ► ^/qS c . 
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-Rmvj [v, /nl] 
-^lv \y, /nl] 



exp 



~S 3 c (a)S 3 (M,a ) 



6cr 2 (M, a ) 



S c (a) 



dS 3 (M,a ) 



1 + -a(M\a )S 3 (M\a ) [5 3 (a) - 35(a)] 



5 c (a)5 3 (M,a )/3 dlogfT 

l d[tr(M|oo)g 3 (M|oo)] 
6 dloger 



1 - i« c (o)S 3 (M, ao) 



5(a) 



5(a) 



(41) 
(42) 



where in the above equations v = S c (a)/a(M,ao), where 
we use the rescaled linear collapse density: S c (a) = 
y/qS c /D(a). In the spherical collapse model, the linearly 
extrapolated density threshold for collapse is 6 C — 1.686. 
Note that, as in the Press-Schechter formalism, we have 
chosen to evaluate the variance and the skewness of the 
density field at the present time ao, and in so doing have 
transferred the time dependence of the theory to the col- 
lapse barrier S c (a) = S c (ao)/D(a). 

Note that these formula are not the same as those orig- 
inally derived by MVJ and LV, but differ by the y/q in 
the definition of the peak height. This was introduced by 
Grossi et al. [20[ in a heuristic way in order to obtain a 
good fit to the FoF mass functions in their numerical sim- 
ulations. They conjecture that it appears for the same 
reason as it appears in the ST mass function [see also 
\Ql\ . The true orginis for such a factor are unclear, since 
for the case haloes identified through a bound spherical 
over-density criterion no q correction is required (l9j | . 

Figure [3] compares the mass function of FoF (b = 0.2) 
dark matter haloes measured from the ensemble of sim- 
ulations at redshifts: z € [1.0,0.5,0.28,0.0] . The left 
and right panels compare the simulation estimates with 
the predictions from the LV and MVJ models, respec- 
tively. The top panels show the absolute mass function 
and bottom panels show the ratio with respect to the 
theory. In these plots we take the Gaussian model for 
the mass function to be that as given by [6l|. Note that 
the errors on the points show the error on the mean, i.e. 
the box-to-box varaince divided by the square root of the 
number of realizations: ct/V iV onscm b ~ cr/3.5. 

A number of important points may be noted. Firstly, 
none of these models fit the data well (26[, however this 
can be mainly attributed to the fact that the ST mass 
function does not fit well the Gaussian simulation data: 
we see a ~10-20% excess in the number of intermediate 
mass haloes M* < M < 100M* and strong suppression 
in the numbers of high mass haloes M > 100M„. Sec- 
ondly, we note that as expected, the model with /nl > 
(/nl < 0) produces an excess (reduction) in the number 
of high mass haloes relative to the Gaussian case. The 
predictions capture these trends. However, as can be 
seen from the lower section of the left figure, the model 
of LV produces the same locus of points for all of the 
models /nl € [±100,0] and for all redshifts. Whereas 
that of MVJ, being almost equally as good, produces a 
slightly different offset for each /nl model and at dif- 
ferent redshifts. This leads us to conclude that in order 



to accurately predict the mass function in non-Gaussian 
models over a wide range of masses and redshifts, one 
simply requires an accurate fit to the Gaussian model, 
combined with the ratio model of LV. 

Figure 3] further emphasizes this point. Here we show 
the fractional mass function for the FoF dark matter 
haloes measured in the simulations at the four expansion 
factors. For each simulation we compute the mass func- 
tion of haloes and take the ratio of this estimate with re- 
spect to the Gaussian model. These results are averaged 
over the 12 realizations and the error are computed on the 
mean. On comparing these estimates with the theoreti- 
cal predictions from the models of LV and MVJ we find 
excellent agreement when the shift parameter 6 C — » y/q5 c 
is used, as advocated in j2£j. There is a small prefer- 
ence to the model of LV, especially at late times for high 
mass haloes in the model with /nl = 100. Coupled with 
the fact that the mass function ratios with respect to the 
ST+LV model are well behaved, we shall hereafter adopt 
the LV model for all our halo model calculations. 

Finally, we note that the ST mass function is preferred 
over other commonly used expressions such as that of 
[ri3 |. since it obeys the important property that when 
integrated over all masses, one recovers the mean matter 
density. 



B. Halo bias 

The halo model also requires us to specify how the 
centers of dark matter haloes of different masses cluster 
with respect to each other. As described by Eq. (|27j) . for 
the Gaussian model, and at first order in the dark mat- 
ter density, halo and matter density perturbations can 
be related thr oug h a scale- independent bias factor b(M). 
Following [13, Hl|, an application of the peak-background 
split approximation enables one to calculate b(M) from 
a given mass function. For the ST mass function the 
Gaussian bias has the form: 

h (\ i i qy2 - 1 i 2p - 6 -M 

bsT{ ">- 1 + ^y + l + ( q ^)p ' V = a(M,z = 0) ' 

(43) 

where the parameters {p, q} are as in Eq. 

As was summarized in [19( . with PNG there are two 
main effects on the bias. Firstly, local non-Gaussianity 
induces a scale-dependent correction factor on extremely 
large scales k < 0.02 /iMpc" 1 [tHHSl- Secondly, as- 
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FIG. 5: Mass dependence of halo bias in models with primordial non-Gaussianity. The left and right panels, show results 
for the halo bias on scales k — {0.005, 0.05} h Mpc -1 , respectively. Top panels: absolute bias vs halo mass. The green, red 
and blue lines correspond to the cases /nl = {0, 100, —100} respectively; dashed lines denote the scale-dependent contribution 
A6 K (fc, /nl); and triple-dot-dash lines denote the scale-independent contribution A&i. Bottom panels: show the ratio of the 
non-Gaussian bias to the Gaussian bias. Solid lines denote the bias ratio for the sum of all components; dashed lines denote the 
same but excluding the scale-independent bias; triple-dot-dash lines denote the same but excluding the scale-dependent bias. 



suming that the peak background split holds, then there 
is also a scale-independent correction to the bias. This 
arises due to the fact that the mass function changes with 
/nl and consequently the typical halo mass changes and 
hence the bias [22,[Z0|- Thus the non- Gaussian bias may 
be written: 

6 NG (£, M, /nl) = b G (M)+Ab K (k, M, / NL )+A6i(M, / nl ) . 

(44) 

The scale-dependent bias term Ab K (k, M, /nl) can be 
written, 



3n m H 2 S c (a) 
2D(a)c 2 T{k)k 2 



A6 K (fc, M, a) = 2/nl [b G (M, a) - 1] 

(45) 

In the original derivation of [l8j . as was pointed out by 
Matarrese and Verde [2j|, a factor of the transfer func- 
tion was missing, and it has been added in the above 
expression. Furthermore, in the original derivation one 
sees that the term in brackets involves just the Gaussian 
bias bc(M, a). However, as we will show in Eq. ([56| . this 
should in fact be the sum of the Gaussian bias plus the 
scale-independent bias correction due to PNG. 



The scale-independent bias correction, A&i, can be 
written 



i gh[%y 

ISOl = 

a av 



(46) 



where R is the fractional correction to the mass function. 
For the case of Lo Verde et al. [60[ the scale-independent 
bias has the form, 



A^ v = - 



1 



6eri?LV 
d{aS 3 ) 
dlnv 



_ dVSa) 
din v 2 



1 



3aS : 



3 \v 



1 



(47) 



It is worth noticing that A6t(/nl) has a sign opposite to 
that of /nl (because the bias decreases when the mass 
function goes up). 

The scale dependence of the above non-Gaussian bias 
model, for haloes with masses above M ~ 10 13 /i _1 M Q , 
was recently tested a gain st a suite of high-resolution TV- 
body simulations by |19l - [21j . For the cross-power spec- 
trum of haloes and matter, it was shown to work at a 
precision of better than 10%. However, the form of the 
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bias has not been investigated for halo masses substan- 
tially lower than M*. In order to make a halo model 
calculation for the mass power spectrum one is required 
to average over haloes of all masses and so it is important 
to understand the behavior for v < 1 as well as for v > 1 . 
Note that whilst the mass function of should only be 
trusted for v > 1 and hence also the bias expansion, since 
we see no deviations at v ~ 1 we shall assume that these 
expansions can be trusted to lower v. 

In Fig.[5J using Eqs (14"T} and (|4"5} , we show how the 
halo bias depends on halo mass, and we extrapolate these 
relations to masses M/M* < 1 (note that in evaluat- 
ing the scale-independent bias we have made the ap- 
proximation d 2 (aS^)/d\nv 2 = 0. Since the bias is also 
scale-dependent, we show the results for spatial scales 
k e {0.005, 0.05} /iMpc -1 . 

On very large scales k ~ 0.005 h Mpc -1 (left panel 
Fig. [5}, we see that the scale-dependent bias correction 
A6 K strongly modulates the Gaussian bias, increasing 
(decreasing) the bias of high mass haloes for /nl > 
(/nl < 0). For low mass haloes with v < 1 this trend re- 
verses. On smaller scales k ~ 0.005 ft, Mpc" 1 (right panel 
Fig. IS}, the scale-dependent bias has been strongly sup- 
pressed and the scale-independent term starts to domi- 
nate. The consequence of this is that high mass haloes 
in models with /nl > (/nl < 0) have a slightly lower 
(higher) bias than in the Gaussian case. Again this trend 
is reversed for haloes with u < 1. 

In evaluating the 2-Halo term we must compute in- 
tegrals of the type as given in Eq. (fT8}. If we consider 
again the large-scale limit U(k\M, /nl) — > 1 as k — > 0, 
then this places the following conditions on the Gaussian 
and non-Gaussian mass functions and bias: 



J dMMn G (M)b G (M) = p ; 
J dMMn NG (M)b NG (M) = p 



(48) 
(49) 

(50) 



These relations arise simply from the fact that 

fi{x,M) = n(M) [l + <S h (x|M)] . 

Taking the non-linear local bias model of [HJ for the 
halo density field, <5 h (x|M) = Y^Lo h(M)5 i (x)/i\, and 
on integrating Eq. (|50} over all haloes, weighted by the 
mass, then we must recover the mass density field: 



dMMn(M) [l + S h (x\M)] 



dMMn(M) 



1 ' 



.(51) 



The integral resulting from the first term in the bracket 
on the right side gives us the normalization condition for 
the mass function / °° dMMn(M) — p, and the integral 
of the second term gives: 



S(x) = - I dMMn(M) V 

P Jo , =n 



bi(M) 



<**(*) ■ (52) 



This can only be true if and only if Eqs (1481) and (|49} are 
true. As a caveat we also note the further condition: 



dMMn(M)bi(M) = (i ^ 1) 



(53) 



Following Eq. (|49} , for PNG this then implies the further 
two conditions that: 



dMMn NG (M) [b G (M) + A6 7 (M|/ NL )] = p (54) 
dMMn NG [M)Ab K {k, M|/nl) = . (55) 



The last condition can only be correct if b G in Eq. (|45 
becomes b G + Abj, i.e. 



Ab K {k, M) = /nl [b G (M) + A6/(M ) - 1] 



Dc 2 T{k)k 2 ' 
(56) 

This correction has not been pointed out in any previ- 
ous work and we expect it to affect the scale-dependence 
of the bias. We note that this should also lead to a 
small asymmetry between the predictions for /nl > 
and those for /nl < 0. 

Returning to the halo model, given Eq. (|55} and the 
fact that as k becomes larger Ab K becomes less important 
(c.f. Fig. |5}, we make the approximation for the dark 
matter, that Ab K = on all scales. For an alternative 
approach see to imposing Eq. ([55} [31]. 



C. Halo density profiles 

1. Gaussian profiles 

The density profiles of dark matter haloes in simula- 
tions evolving from Gaussian initial conditions has been 
studied in great detail. A reasonably good approximation 
for the spherically averaged density profile is the Navarro, 
Frenk & White model !71]. This can be written: 



pNFw(r|M) = p5 c {M) 



r 
ro 



r 
ro 



(57) 



where the two parameters are the scale radius ro and the 
characteristic density 5 C . Note that if we define the halo 
mass to be M v i r = 47rr^ ir 200p/3, then owing to mass 
conservation there is only one free parameter: the con- 
centration parameter c(M) = r v ; r /ro, and we have 



S C (M) 



200c 3 /3 



log(l + c) - c/(l + c) 



(58) 



The parameter c(M) can be obtained from the original 
model of NFW, but instead we prefer to use the model 
of Bullock et al. [72j]. Note that we correct c(M) for 
the fact that the definitions of the virial radius in the 
Bullock et al. model for the concentration parameter and 
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FIG. 6: Ensemble average density profiles of the dark matter haloes in the Gaussian (/nl = 0) simulations as a function of 
radius. The panels show haloes with bin centered masses log 10 M//i _1 M ] € {13.62, 14.22, 14.92, 15.42}. In all panels: points 
with errors show estimates from the simulations; red solid lines give the results for the NFW model; blue solid lines give results 
for NFW model convolved with a Gaussian filter of radius 2.5 times the softening scale. The dotted lines either side of the main 
prediction, give the predictions if the halo mass associated with the upper and lower edges of the mass bin are used. Vertical 
dashed lines denote the softening length and vertical dot-dashed curves denote r v i r for the central mass of the bin. Note, the 
density profiles are constructed from only those partcicles that are members of the FoF groups. 
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FIG. 7: Ensemble average of the ratio of the density profiles in the non-Gaussian model with the profiles in the non-Gaussian 
models as a function of radius. The red and blue colors corresponding to /nl = 100 and /nl = —100, respectively. The points 
with errors denote estimates from the simulations; solid lines denote the log-linear profile ratio model as given by Eq. (|61|) . and 
the vertical lines are as in Fig. [6] 
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TABLE I: Expected number of haloes in the mass bins from 
which density profile averages are calculated, per simulation. 
Columns are: (1) number of mass bin; (2)-(3) lower and upper 
edges of the mass bin; (4)-(6) number of haloes in bin. 



Bin 

# 


logio M i 

[h- x M & ] 


logic M 2 


# haloes 

/nl = 


# haloes 
/nl = 100 


# haloes 
/nl = -100 


1 


13.169 


13.469 


617526.6 


613869.9 


620691.9 


2 


13.469 


13.769 


311273.1 


309856.2 


312781.2 


3 


13.769 


14.069 


144713.0 


144526.5 


144994.3 


4 


14.069 


14.369 


60805.3 


61050.2 


60481.0 


5 


14.369 


14.669 


21744.9 


22108.4 


21371.4 


6 


14.669 


14.969 


5843.5 


6080.4 


5608.2 


7 


14.969 


15.269 


1008.0 


1092.5 


922.5 


8 


15.269 


15.369 


86.7 


101.1 


72.7 



the Sheth & Tormen mass definition used for the mass 
function are different [for details as to how to do this see 
H§| . Over the past decade a number of alternative models 
for halo profiles have emerged. Owing to the relatively 
low resolution of our haloes we believe that the original 
model of NFW will be of sufficient accuracy to describe 
our haloes. 

Figure [5] shows the ensemble average density profiles 
of dark matter haloes in the simulations. The haloes 
were separated into a set of mass bins of equal loga- 
rithmic width Alog 10 M = 0.3, and with the minimum 
halo mass from which a profile can be estimated being 
taken as 50 particles (M ~ 1.5 x 1O 13 /i _1 A/ ). Note 
that, whilst the number of particles is relatively small 
for the lowest mass haloes used in the profile estimation, 
as Table U shows, we are averaging over large numbers 
of haloes and multiple simulations per mass bin. The 
figure shows the results for a sub-sample of four of the 
mass bins, and these correspond to bins {2,4,6,7} = 
{13.62, 14.22, 14.92, 15.12} [log 10 M/hr x M Q ] in TableU 

We estimate the density profile for each individual 
halo by taking only the particles that are in the FoF 
halo (b = 0.2). We compute the halo center of mass 
and the radial distance of each particle from this cen- 
ter. The particles are then binned into equal loga- 
rithmic radial bins of thickness, A log 10 r[ /i _1 Mpc] 
(I, I Our estimate for the profile is then given by 
p{n,M) = mpNijlr^Mjl/Vehieiiir), where N^n, Mj] 
is the number of particles in the ith radial bin for 
haloes in mass bin Mj and the shell volume is V s h c \\ — 
4n [(n + Ar,/2) 3 - (r< - Ar,/2) 3 ] /3. 

In Fig. [5] we see that the model and the data do not 
agree at the inner and outer parts of the profile. The data 
are significantly flatter in the inner radius than the NFW 
model would suggest. However, the softening length for 
the simulations was r ~ 0.04 h Mpc, and so at this scale 
we expect the core to be effectively of constant density. 
The vertical dash line in each plot shows the softening 
length. The virial radius taken plotted for the 

mean mass in the bin and this is the vertical dot-dash line 
in each panel. There is reasonable agreement between the 
scale at which the FoF haloes are truncated and r v ; r . 



We considerably improve the agreement between the 
NFW model and the data by convolving the theoretical 
profiles with a Gaussian filter function of radius 2.5 times 
the softening length, i.e. 



PNFw(r\M) 



d 3 k 

(2tt) 3 



MU NF w{k\M)W(k)j (kr) , 



(59) 

where W(k) = cxp [— (2.5/ so ftfc) 2 /2] and where 
U(r\M) = p(r\M)/M is the mass normalized profile. We 
truncate the profile at the virial radius and for the NFW 
model the Fourier transform can be written as (481 ] : 



/(c)C/ NFW (fc|M) = - 
+ cos [kr a(l 



sin(/ccTo) 



fcr (l + c) 
c)]{Ci[kro(l + c)]-Ci[kr ]} 

sin [fcr„(l + c)] {S t [fcr (l + c)] - S,:[fcr ]}(60) 



where /(c) = [log (1 -I- c) — c/ (1 + c)] and where Si and 
Ci are the standard sine and cosine integrals. This makes 
the agreement between the low-mass halo samples and 
the data almost perfect, however the higher mass halo 
samples would require larger smoothing radii to explain 
the difference. One might motivate this by the fact that 
for the highest mass bins these haloes are just forming 
and as such are more likely to have complex structure, 
and so the center of mass may not be a good proxy for the 
halo center. A better choice may be the point of deepest 
potential. We shall not pursue this matter further, but 
note that for the highest mass haloes in our simulations 
the profiles appear to be flatter inside r ~ 0.5 h Mpc. 



2. Primordial non-Gaussianity and density profiles 

The impact of PNG on the density profiles of dark 
matter haloes has not been investigated previously. We 
therefore make a short exploration. 

Fig.[7]shows the ratio of the density profiles in the non- 
Gaussian simulations with those from the Gaussian ones. 
We estimate this quantity by dividing the mean halo pro- 
file for each mass bin for each /nl simulation with the 
corresponding one from the Gaussian runs. The plotted 
points with errors are then the average and 1-a errors 
from the 12 simulations. In the figure we clearly see that 
there is an effect of PNG on the profiles. For the case 
of /nl > 0, we find that the profiles are denser in the 
inner regions of the halo; the converse is true for the case 
/nl < 0. For /nl = ±100, the strength of the effect 
depends on the halo mass considered: for cluster, group 
and small group mass haloes we find effects of the order 
{< ±4.5%, ± < ±3.5%, < ±2.5%}. 

Since the power spectrum in the halo model depends 
on both the profile and the square of the profile, these 
effects are important to characterize for accurate cluster- 
ing predictions on small scales. We therefore attempt to 
model this in a simple way. From Fig. [7] it can be seen 
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FIG. 8: Best fit parameters for the density profile ratio 
-Rp(/nl, M), as defined in Eq. (|61[) . as a function of the FoF 
halo mass, measured at z — 0. Top panel: the zero-crossing 
parameter rx- Bottom panel: the slope parameter m. The 
solid red lines denote the results for /nl = 100 and the blue 
dashed lines the results for /nl = —100. 

that the profile ratio as a function of log-radius appears 
to be almost a straight line. We therefore fit a log-linear 
model to the measured profile ratios. Explicitly the ratio 
model has the form: 

R p (r, M, /nl) - 1 = m log 10 [r/r x ] (61) 

where m is the change in the slope and rx is the zero 
crossing scale, both of which depend on halo mass and 
/nl- We fit for the parameters {m, rx} over the range of 
radii (Z so ft < r < r v iri)) and the resultant model fits are 
shown in Fig. [7] as the solid lines. The best fit parameters 
as a function of halo mass are shown in Fig. [31 where we 
see that, for all but the lowest mass bin in the /nl = 100 
model, the values of rx increase with increasing mass and 
that these are almost identical for both the positive and 
negative /nl models. On the other hand, the values for 
the slope m monotonically decrease/increase for /nl pos- 
itive/negative. Modulo the sign, these values are similar 
for both /nl models. 

To use this correction ratio in the halo model, we spline 
fit {rx, m} as a function of mass, with the exception of 
the lowest mass bin. We assume that the asymptotic 
limit for low masses is R p — > 1 as M — > 0, and enforce 
this by adding {r x ,m} = {0,0} for M = 10 6 h~ 1 M Q as 
an extra data point in the spline fitting. 

What is actually required for the calculation of the 
power spectrum is the mass normalized Fourier trans- 
forms of the density profiles, and we show this in Fig. [H] 
We calculate the Fourier transform of the profile as in 
Eq. (|59p . Unfortunately, owing to the radial dependent 
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FIG. 9: Dependence of the Fourier transform of the mass 
normalized density profile on scale for several halo masses. 
Top panel: absolute value of the Fourier transform of the pro- 
file (U(k\M) — p(k\M)/M) as a function of wavenumber in 
units of inverse virial radius. The solid green, dashed red and 
dot-dashed blue lines show results for /nl £ [0, +100, —100]. 
The curves from top to bottom show results for haloes with 
masses M G [5.0 x 10 13 , 1.0 x 10 14 , 5.0 x 10 14 , 1.0 x 10 15 ,5.0x 
lO 15 ]/i _1 M0. Bottom panel: ratio of the non-Gaussian to 
Gaussian profiles in Fourier space. The higher mass haloes 
show stronger amplification/suppression for /nl = 100/— 100. 



correction factor Eq. (|61|) . there is no analytic solution 
to this integral and so we compute this numerically. We 
do however make the following improvement to compu- 
tational speed: the profile can be rewritten 



U(k\M, / NL ) = U(k\M)R p [k, M, /nl] (62) 



where U(k\M) is given by Eq . (f60|) . We then gener- 
ate a bicubic spline [see |73|, for details] fit to the ra- 
tio Rp as a function of kr v - u - and halo mass M. As 
a final note, the above expression is not normalized 
correctly in that U(k\M, /nl) does not approach unity 
on large scales. We are however free to renormalize 
U{k\M, /nl) and this can be done through the opera- 
tion: U(k\M, / NL ) - U(k\M, f NL )/U(k = 0|M, / NL ). 
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FIG. 10: Comparison of the matter power spectrum in the halo model and measurements from an ensemble of numerical 
simulations (see £11111 for details). Top panel: Absolute power. Points with error bars show estimates from the simulations. 
The halo model predictions are: red dash line denotes Pih; blue dot-dash line denotes P2H; magenta dotted line represents 
Pz£ c \ the black solid line shows the total halo model prediction including subtraction of the halo exclusion term; the green 
triple-dot-dash curve is the same, but neglecting the halo exclusion term. Bottom panel: ratio of the matter power spectra with 
predictions from halof it 56]. Points and line styles are as above. 



VI. RESULTS 

We now put together all of the components of the halo 
model and make predictions for the nonlinear matter 
power spectrum and the matter correlation function in 
the Gaussian case and then for the models with local 
PNG. A practical note on evaluating the 2-Halo term 
for the mass distribution: we are required to compute 
integrals over an infinite domain in halo mass. This is 
computationally challenging, and so instead we make the 
following approximation: 



P. 



2H 



Pi 



NL 



where p = dMn{M)Mb(M). Provided M x and M 2 
are sufficiently small and large halo masses, then the 
above integral approaches the exact answer of infinite 



M-2 



P J Mi 



dMn(M)Mb(M)U (k\M) 



(63) 



limits over a restricted range of k. We set {Mi,M2} 
{1.0 x 10 6 ,1.0 x lO^j^M©. 



A. Estimating the power spectrum 



The density Fourier modes were estimated using the 
conventional Fast Fourier Transform (FFT) method: the 
dark matter particles were assigned to a regular cubi- 
cal grid using the 'Cloud- In-Cell' (CIC) scheme [74[, and 
throughout we used N g = 1024 3 sized Fourier meshes. 
The FFT of the gridded density field was then computed 
using the FFTW routines [7|| . Each resulting Fourier mode 
was corrected for the convolution with the mesh by di- 
viding out the Fourier transform of the mass-assignment 
window function. For the CIC algorithm this corresponds 
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k [Mpc/h] k [Mpc/h] 

FIG. 11: Comparison of the matter power spectrum in models with Gaussian and non-Gaussian initial density fluctuations at 
redshifts z = 1.0 (left sub-fig.) and z = 0.0 (right sub-fig.). Top panels: Absolute power. Points with error bars show results for 
the simulations and the colors green, red and blue denote the models /nl = {0, +100, —100}. The lines represent halo model 
predictions: dash lines denote Pih; dot-dash lines denote P2H; dotted lines denote P^i c ; the solid line represents the total halo 
model prediction including subtraction of the halo exclusion term. Bottom panels: ratio of the matter power spectra in the 
/nl = +100 and —100 models with respect to the Gaussian (/nl = 0) results. Points and line styles are as above. 



to the following operation: 



sin [7r/cj/2fc 



NyJ 



[7rfci/2fc Ny ] 



(64) 

where sub-script d and g denote discrete and grid quan- 
tities, and where k^y = irN g / L is the Nyquist frequency. 

The discrete power spectra on scale fc; are then esti- 
mated by performing the following sums, 



(65) 



1=1 



where Nk is the number of Fourier modes in a spherical 
shell in fc-space of thickness Afc. 



B. Matter power spectrum: Gaussian case 

Figure [TU] presents the nonlinear matter power spec- 
tra measured in the simulations for the Gaussian case 



at redshifts z = 1.0 and z — 0.0, left and right panels 
respectively. The top sections of the panels show the ab- 
solute power and the lower ones the ratio with respect to 
the predictions from halof it [5^. The figure also shows 
a term by term breakdown of the halo model predictions. 
It can clearly be seen that the sum Pih + Pzh without 
subtracting the exclusion term over shoots the measured 
nonlinear power by a factor of 2-3 (green triple dot dash 
line). The subtraction of the term due to halo exclusion 
Pih + P-m — Pffl° significantly improves the predictions 
of the model on all scales and at both of the redshifts 
considered. We see that at z = 0.0 the total halo model 
result is within a few percent of the measured power from 
the numerical simulations on all scales measured. At 
z = 1.0 the same statement is true except on small scales, 
k > 0.5/iMpc -1 where there appears to be a significant 
disagreement of the order 20%, this is partly due to the 
failure of halof it on these scales (which is the input 
nonlinear power spectrum for the 2-Halo term). Improv- 
ing the modeling of the halo-halo center clustering will 
most likely solve this issue. This discrepancy may also 
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FIG. 12: Matter correlation function as a function of spatial scale measured in the Gaussian simulations at redshifts z — 1.0 (left 
sub-fig.) and z — 0.0 (right sub-fig.). Top panels: absolute correlation function. Points with errors denote the measurements 
from the simulations and the lines denote the Halo Model predictions: dash lines denote £ih; dot-dash lines denote ^2h; dotted 
lines denote ^2H C ; solid lines represent the total halo model prediction including subtraction of the halo exclusion term; triple 
dot-dash lines show the same but without subtracting the exclusion term. Bottom panel: ratio of the measurements and halo 
model predictions with respect to the nonlinear correlation function from halof it. 



be alleviated by the fact that the Nyquist frequency for 
the power spectra is k^ y = 7rl024/1600 ~ 2.0 h Mpc -1 , 
and that one should only really trust the results up to 
fc Ny /2 - 1.0/iMpc" 1 . 



On large scales we see that the Pm term 
(dashed red line) asymptotes to a constant value ~ 
{40, 400} h~ 3 Mpc 3 at z = 1 and z = 0.0, respectively. 
This is significantly larger than the expected amplitude 
due to the shot-noise of the particles 1/n = V/N ~ 
3.8 h~ 3 Mpc 3 . The term (magenta dotted line) ef- 

fectively kills the excess shot-noise of the 1-Halo term. 
This is entirely consistent with the theoretical expecta- 
tions in Appendix [Bj The subtraction of the term due 
to halo exclusion therefore is an essential correction to 
make in order to make realistic predictions in the halo 
model. 



C. Matter power spectrum: PNG case 

Figure QT] shows the effects of PNG on the matter 
power spectrum at redshifts z — 1 and z = 0, left 
and right sub-figures respectively. The measurements 
with errors are the ensemble average power obtained 
from the 12 realizations per /nl- The predictions from 
the halo model calculation are also plotted and we see 
good agreement. However, the difference between the 
/nl = {0, +100, —100} models is very small and so we 
take the ratio of the PNG models with respect to the 
Gaussian case, and this is what is plotted in the lower 
sections of the figures. This clearly shows, as was seen 
earlier in Fig. [2] that the changes are of the order ~ 3.5% 
at z = 1.0 and of the order ~ 2.5% at z = 0.0. The halo 
model predictions for this ratio are in excellent agreement 
with the measured ratios. 

The strength of the effect appears to peak around 
k ~ 1, and then declines at higher k. This can be un- 
derstood in the following way: as one goes to higher k 
the 1-Halo term comes to dominate. The integrand for 
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FIG. 13: Comparison of the matter correlation function in models with Gaussian and non-Gaussian initial density fluctuations 
at redshifts z — 1.0 (left sub-fig.) and z = 0.0 (right sub-fig.) Top panels: absolute correlation function. Points with errors show 
results for the simulations and the colours green, red and blue denote the models /nl = {0, +100, —100}. The lines represent 
halo model predictions: dash lines denote £ih; dot-dash lines denote ^2h; dotted lines denote Cff-T; the solid line represents the 
total halo model prediction including subtraction of the halo exclusion term. Bottom panels: ratio of the matter correlation 
functions with /nl = +100 and —100 with respect to the Gaussian (/nl = 0) results. Points and line styles are as above. 



this is n(M)M 2 U(k\M), and for cluster mass scales it 
peaks around k ~ 1 and then decays strongly, whereas 
for group and galaxy mass scales it peaks at much higher 
k. As we have seen earlier in S]V] the effect of PNG on 
the mass function and density profiles is strongest for 
the highest mass haloes. This then leads us to expect 
that the effect of PNG on the nonlinear power spectrum 
should peak around k ~ 1. The dashed lines in Fig.fTTI 
show the 1-Halo term and so give confirmation of this 
logic. 



D. Matter correlation function: Gaussian case 

We estimate the matter correlation function for the 
2 = 1 and z — snapshots in each simulation to an 
accuracy of order 3-5% using our fast and exact corre- 
lation code DualTreeTwoPoint, which is based upon the 
fcD-tree data structure, and the code is parallelized using 
MPI calls. Thus on averaging over the 12 simulations we 
expect results that are accurate to 5%/Vl2- 1 < 2%. 



Figure presents the ensemble average estimate of 
the matter correlation function in the Gaussian models 
over three decades in spatial scale at redshifts z — 1.0 and 
z = 0.0, left and right panels respectively. The figure also 
shows the halo model predictions appear in remarkably 
good agreement with the simulation data. The exact 
deviations are hard to quantify on the log-scale and so we 
take the ratio of the theory and simulation measurements 
with respect to the halof it model correlation function. 

We now see that the halo model predictions are better 
than 10% over the entire range of scales and redshifts con- 
sidered. The predictions are somewhat worse at the 2- 
to 1-Halo cross-over scale (i.e. r £ [2, 10] /i _1 Mpc), also 
on the very largest of scales around the BAO feature and 
on the smallest scales r < 0.2/i _1 Mpc. We emphasize 
that none of the halo model parameters were tuned to fit 
the clustering statistics directly. 

In the figure we also show the result for the halo model 
calculation if no exclusion correction is made, and we 
see that predictions significantly overshoot the measure- 
ments by factors of a few on small scales, especially at 
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low redshift. The figure also shows that the exclusion cor- 
rection essentially kills the contribution from the 2-Halo 
term to the correlation function on small scales. Further- 
more this correction also kills some of the contribution 
of the 1-Halo term to the correlation function on scales 
larger than r ~ 2/i _1 Mpc. 



E. Matter correlation function: PNG case 

In Figure [13] we present the ensemble average estimate 
of the matter correlation function in the models evolv- 
ing from PNG initial conditions at redshifts z = 1.0 and 
z = 0.0, left and right panels respectively. As for the 
correlation function in the Gaussian case the halo model 
predictions with exclusion provide a remarkably good de- 
scription of the clustering. The differences are not clearly 
visible on the log-scale, and so we take the ratio of the 
PNG models with the Gaussian case. Note that we con- 
struct this ratio for each simulation and compute the 
ensemble average. 

The lower panels of the figures show that, as in the 
case for the matter power spectrum, there is a significant 
signal of /nl on the small-scale correlation function. At 
z = 0.0 the signal is of the order ~ 2.5% and affects the 
clustering on scales r < 5h~ 1 Mpc. At higher redshift, 
z = 1, we clearly see the same general trend but the 
relative difference in the signal is a little larger ~ 3 — 4%. 
However, the measurements appear slightly noisier. 

Once again the halo model predictions with exclusion 
provide an excellent description of the ratio, being accu- 
rate to 1% precision. 



VII. CONCLUSIONS 

In this paper we have developed the halo model of 
large-scale structure for application to computing mat- 
ter clustering statistics in models with PNG initial con- 
ditions. In particular, we considered the case of local 
quadratic corrections to the primordial potential, char- 
acterized by the parameter /nl, see fJTT] for details. 

In ijllll we provided details of the large ensemble of 
A-body simulations that we employed in this study. 

In jjlVI we explored standard nonlinear perturbation 
theory techniques to predict the matter power spec- 
trum in models with PNG. It was demonstrated that the 
next-to-leading-order correction to the power spectrum 
worked very well up to scales k < 0.2 h Mpc -1 for the ra- 
tio with respect to the Gaussian case, but that on smaller 
scales the PT failed to reproduce the simulation results. 
Further, the absolute power was only well reproduced by 
the PT on scales < 0.1 hMpc^ 1 . 

In mVl we reviewed the halo model and gave a calcu- 
lation of halo exclusion. In Appendix [B] we also showed 
theoretically that exclusion can help resolve the problem 
of the large-scale excess power in the halo model. 



In SjV] we performed a numerical study of the ingredi- 
ents of the halo model, in the context of PNG. We stud- 
ied the halo mass function; the halo bias and the density 
profiles. For the mass function ratios, we confirmed that 
the existing models of Lo Verde et al. [60( and Matarrese 
et al. fl3j were in excellent agreement with the simula- 
tions, once the peak height was rescaled [2(3]. We found 
that the Gaussian mass function of Sheth and Tormen 
[6l| was a poor description of the simulation data. 

In |V]we examined the halo bias in the context of PNG 
and showed that, if the halo model is correct, then there 
must be a small asymetry in the scale-dependence of the 
bias in non-Gaussian models on large scales. 

In SjV]we explored the density profiles of dark matter 
haloes in the context of PNG. This has not been under- 
taken before and we showed that halo profiles become 
more (less) dense in the presence of /nl > (/nl < 0). 
We found, for /nl = ±100, that cluster, group and 
small group mass haloes were modified at the level of 
{< ±4.5%, ± < ±3.5%, < ±2.5%}. We modeled these 
effects with a simple log-linear model. 

In i j VII we presented the two-point clustering statistics 
of the matter in the simulations and in the halo model. 
We showed that including halo exclusion in the model 
was important in order to produce accurate results. We 
also demonstrated that on large scales the halo exclusion 
term almost exactly canceled the excess large-scale power 
arising in the 1-Halo term. This appears to solve a long 
standing technical problem for the halo model. 

For the case of the absolute power spectra, whilst the 
theory and the measurements differed by up to 10%, the 
ratio of the non-Gaussian to Gaussian predictions and 
measurements differed by of the order ~ 1%. For the 
case of the correlation functions, the absolute predic- 
tions on small scales were very good and the ratio of the 
non-Gaussian to Gaussian predictions and measurements 
were also accurate to ~ 1% 

We conclude that the modeling that we have developed 
in this paper, will be good enough for predicting the abso- 
lute value of matter clustering statistics to within < 10% 
and their ratios to < 1%. We anticipate that this will 
be useful for constraining /nl from measurements of the 
nonlinear shear correlation function in future weak lens- 
ing surveys such as EUCLID [H and LSST j3|. 
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[78] At first glance the nonlinear (in /nl) corrections to the 
power and bispectrum in Eqs © and ©, appear to be 
infra-red divergent quantities. This can be seen most sim- 
ply by taking the limit k — > 0, and for the power spec- 
trum this leads to an integral of the form J °° dqq 2 P 2 (q). 



Taking -P</>(q) oc k" 4 , then the above integral becomes 
/ °° £% 2(n ~ 4)+2 oc k a \™, where a = 2(n - 4) + 3. Taking 
n — 1, then we get a — — 3, and we now see that if we 
substitute the lower limit then the expression diverges. 
However, as was pointed out by [26|, if one believes that 
inflation was responsible for the non-Gaussianity, then 
one should invoke the cut-off scales which limit the phys- 
ical processes: on small scales the fluctuations are cut-off 
by the reheating scale during inflation, and on the large 
scales, by the horizon size today. In this case and for re- 
alistic values of /nl, then the higher order corrections 
should be subdominant at the present time, and so can 
be neglected [H,[ll. 
[79] Note that we assume that haloes, which are found with 
FoF halo finders in numerical simulations, have a mass 
Mpop and that this corresponds to Mf f = M200 = 
4/37rr|oo200p, i.e. the mean overdensity inside a sphere 
of radius r2oo is 200 times the mean background den- 
sity. The only issue that one must be careful of when 
dealing with halo phenomenology is that the same mass 
convention is adopted for all semi-analytic models: mass 
function, profiles and bias etc. Where there are clashes 
of definition, we are careful to obtain mappings between 
these schemes to ensure consistency following the method 
of Smith and Watts 
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Appendix A: Some results concerning the primordial Skewness 



In this appendix we give some relations concerning the skewness of the density field. The ensemble average of the 
cube of the density field using Fourier space representation can be written: 

d 3 k! d 3 k 2 d 3 k 3 



<<&(*)) = 



(2tt) 3 (2tt) 3 (2tt) : 



-W{k u M)W(k 2 ,M)W(k 3 , M) (5(ki)5(k 2 )5(k 3 )) 



— i(ki+k2+k3)-x 



(Al) 



For a Gaussian field the average of the three Fourier modes gives zero, but for the local model of non-Gaussianity the 
product is related to the primordial potential bispectrum. Using Eqns ([5]), and ©, we arrive at, 



(S 3 m(x)) = 



d 3 ki d 3 k 2 d 3 k 3 



(2tt) 3 (2tt) 3 (2tt) 3 
x ($(k 1 )$(k 2 )$(k 3 )) exp[-i(k! + k 2 + ka) • x 
d 3 ki d 3 k 2 d 3 k 3 



W(k u M)W(k 2 , M)W(k 3 , M)a(k 1 )a{k 2 )a(k 3 ) 

i)exp[-i(k a +k 2 + k 3 ) • x] 

( .,_ ) , ( ,^ ) , ,,, W(ki,M)W(k 2 ,M)W(k 3 ,M)a(ki)a(k 2 )a(k 3 ) 

xB, I ,(k 1 ,k 2 ,k 3 )(27r) 3 (5(k 1 +k 2 + k 3 ) exp[-i(kx + k 2 + k 3 ) • x] 
d 3 k d 3 k d 3 k 

^j^j^W(ki,M)W(k 2 ,M)W(k 3 ,M)a(k 1 )a(k 2 )a(k 3 ) 

x2/ NL [P <A (fc 1 )P < ^(fc 2 )+2cyc] (2 7 r) 3 (5(k 1 + k 2 +k 3 )exp[-t(k 1 + k 2 + k 3 ) • x] . (A2) 

Computing the delta function integrals separately for each of the three terms in square brackets and defining 
k 2 = k\ + k\ + 2kik 2 fi, we then find that the skewness (S3 = (<5 3 ) / (<5 2 ) ) has the form: 



S 3 a 4 = 3/, 



XL 



d 3 ki d 3 k 2 



M(7 1 -|.)P,,(/,- l iM|/ 1 -„.i/',,(/,-. 2 ) / dfia(k 3 )[W(ki,M)W{k 2 ,M)W(k 3 ,M)] . (A3) 



(2tt) 3 (2tt) 3 

We are also interested in computing the derivatives of the skewness with respect to the mass variance a up to 
second order: 

d(S 3 a 4 ) 1 dlogM d(S 3 <r 4 ) 



3 d log a dlogR 

1 d 2 logM d(S' 3 CT 4 ) 1 fdlogM\ 2 d 2 (S 3 a i ) 



d log a 
d 2 (S 3 a) 

d(logcr) 2 3d(loger) 2 dlogR ' 9\d\ogcrJ dlogR 2 ' 
where we used the fact that d log M/d log R = 3. The required auxiliary functions are: 



(A4) 
(A5) 
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NL 
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NL 
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dlogff 2 (P) 



dlogM 



(A6) 

(A7) 
(A8) 

(A9) 



dlogAf 2 9ct 2 (P) J (2tt) 3 ^ 

where we have introduced the notation W l = W{k % R), W[ = dW {hR) / dlog(hR) and W[' = d 2 W {hR) / dlog^R) 2 . 
For a real-space top hat filter function we have: 

3 

Wth(v) = — [sin y -y cosy] ; (A10) 



y 

3 

Wtk(v) = ^ [(y 2 ~ 3)sin?/ + 3ycos2/] ; 



[(9 - Ay 2 ) sin y + y (y 2 - 9) cos y] 



(All) 
(A12) 
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Finally, some relationships that will be of use for calculating the halo mass function and also the scale-independent 
contribution to the halo bias are: 



dlog(aS 3 ) 
dlogv 

d 2 log(aS 3 ) 
d log v 2 



1 



<J A S 3 



3a 4 S, 



d log a 
d^Ss) 



d 2 (a 4 S 3 )- 



d log a d log a 



(A13) 
(A14) 



Appendix B: A possible solution to the problem of excess power on large scales in the halo model 



Several authors have pointed out that the halo model fails to reproduce the correct clustering statistics on large 
scales p8l . l5ol l59l [76| • This arises because the 1-Halo term approaches a constant value significantly in excess of the 
shot noise for dark matter particles in simulations: e.g. 



i r°° i 

lim P m (k) = ttt / dMn(M)M > — , 
fe^o p z J n s 



(Bl) 



where n$ is the number density of dark matter particles. A more recent criticism of the halo model, in the context of 
renormalized perturbation theory (RPT), was given by (59j . We now show how halo exclusion may resolve this issue. 
Consider the RPT formulation of the nonlinear matter power spectrum in perturbation theory [59| : 



iWfc) = G 2 (k)P Un (k) + P MC (k) 



(B2) 



where G(k) is the nonlinear density propagator, which informs us of how a density mode decorrelates from its initial 
state; and the term Pmg informs us of the power gained by a single mode from coupling with all other modes. (59j 
make the analogy P 2 h — > G 2 (k)Pun and Pm — > Pmc(&)' They then measure the quantity P ~~ G 2 (k)Pu n from 
simulations and compare the result with the theoretical predictions for Pmc ^ k 4 and Pm oc const. The simulations 
show a fc 4 slope, and hence they conclude that the halo model fails on large scales. 

Let us now reexamine this issue in the context of halo exclusion. Firstly, consider again the halo center clustering 
of haloes of masses Mi and M 2 , if we treat this as in (sijj and expand the halo density as a local function of the 
nonlinear dark matter density then we have: 



P 



hh 



6 1 (M 1 )6 1 (M 2 )P NL (fc|#) + 0(62,63,. 



(B3) 



Consider now the very large-scale limit of the halo model, then we have U — > 1 as k — > 0. Recalling the conditions 
Eq. (gS]), Eq. (gSJ) and Eq. ([52]). the halo model then reduces to, 

Pnu(k) = PMk\R) - P 2 °H C (k) + Pm(k) . (B4) 

If we now insert the RPT expansion for the large-scale smoothed nonlinear matter power spectrum, then we have: 

Pau(k) = \W(k\R)\ 2 [G 2 (k)P Un (k) + P MC (k)] - P^ c (k) + P m (k) . (B5) 

On rearranging the above equation, one hnds 

P HM (fc) - \W(k\R)\ 2 G 2 (k)P Un (k) = \W(k\R)\ 2 P MC (k) - P^ c (k) + P m (k) . (B6) 

In the limit k ^ 0, W(k\R) — > 1. The right hand side of the above equation can be consistent with the simulation 
results of [59j . if and only if Pm — P^ — > e, where e is a sufficiently small quantity that the right side of the above 
equation decays as oc fc 4 for all observable scales of interest. Using Eq. (f34|) we may explicitly write e as: 



dMr 



Y[ I -_ J dM in {Mi)M^ V(n + r 2 ) [1 + 6(Mi)6(M 2 )e N L(n + r 2 )] . (B7) 



On taking the relation between halo volume and mass to be M v ; r = 47rr l y lr Aj(?/3, where A is an overdensity threshold 
that defines the halo today, e.g. A ~ 200, then we may expand the exclusion volume in powers of the halo mass, 



V{r 1 +r 2 ) 



1 

^A 



Mi + 3M 2/3 M?, /3 



3All /3 M 2/3 



My 



(B8) 
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On inserting this into the second term on the left side of Eq. (|B7I) we find, 



f 



p 3 A 



dM 1 n(M 1 )M 1 I dM 2 n(M 2 )M 2 1M X + 6M* /3 M 2 /3 [l + b(M 1 )b(M 2 )^ N L(ri + r 2 )] 



(B9) 



Let us consider two possible limits of the above expression: 
• £nl ^ 1 ; then we have 



p 3 A 
2 



dM 1 n(M 1 )M 1 / dM 2 n{M 2 )M 2 [2M1 + 6Mi /3 M 2 1/a 



dM x n{M x )Ml + 



p 3 A 



dMm(Mi)M 1 5/3 / dM 2 n(M 2 )M, 



.4/3 



• £nl ^ 1 ; then we have 



p 3 A 



dMin(Mi)Mi&(Afi) / dM 2 n(M 2 )M 2 b(M 2 ) + 6M 2/3 M 2 /3 ] £nl(>i +r 2 ) 



(BIO) 



(Bll) 



To get an estimate of this, suppose that the first term in Eq. (jBllI) dominates over the second and that 
max [^nl^i + r 2 1 -R) ] ~ /A, then we would have 



dMn(M)M 2 [1 - 2fb(M)] 



(B12) 



We thus see that the correction can be very close to the order of the resulting large-scale power for the 1-Halo term, 
and depends strictly on the quantity 1 — 2/6(M). We note that the case / = 1 would mean £ ~ A, and that this 
would be unrealistically large. This may be argued in the following way: A is the volume average density at the viral 
radius, whereas bi(Mi)bi(M 2 )£ iNL (ri + r 2 \r) is the volume average correlation function of halo centres smoothed on 
the scale R, for separations n + r 2 , and excluding the points within the same halo, thus we conclude that / < 1. In 
fact we would need max [£ NL (ri + t^I-R)] < A/26(M). 

As we show in Fig.[TUl exact evaluation of the 1-Halo and the 2-Halo exclusion terms on large scales produces two 
quantities that almost exactly cancel. We therefore forward halo exclusion as a risible solution to the problem of 
excess large-scale power in the halo model. 

As this paper was being submitted to the arXiv, a paper discussing an alternative solution to resolving the break 
down of the Halo Model on large scales was proposed by [77 \. 



